Estimated Costs and Outcomes Associated With Use and Nonuse of Medications for Opioid Use Disorder During Incarceration and at Release in Massachusetts

Key Points Question Is provision of medications for opioid use disorder (MOUD) during incarceration associated with fewer overdose deaths? Findings This economic evaluation of a model of the natural history of OUD in Massachusetts found that a strategy offering buprenorphine, methadone, and naltrexone during incarceration was associated with 192 fewer overdose deaths (a 1.8% reduction) and was less costly than a naltrexone-only strategy averting 95 overdose deaths (a 0.9% reduction). The 3-MOUD strategy was also cost-effective at $7252 per quality-adjusted life-year gained. Meaning These findings suggest that offering 3 MOUDs during incarceration is a life-saving, cost-effective intervention.

The user can create either an open or a closed cohort simulation. In an open cohort simulation, new population "arrives" to the simulation in every time step such that the total population in the model always reflects the size of the total population with OUD living in that jurisdiction. The arrival rate represents both the development of new opioid use disorder and migration into the state among those with existing opioid use. In a closed cohort, no cohort members enter the simulation and the size of the population in the simulation dwindles over time as cohort members die.
The core simulation (Supplemental Figure 1) of the RESPOND model involves the simulation of the natural history of opioid use disorder as a relapsing and remitting disease over a lifetime. RESPOND simulates OUD as a series of transitions between four health states of opioid use: 1) active, non-injection, 2) non-active, non-injection, 3) active injection, and 4) non-active, injection opioid use. In each time-step of the simulation, population fractions move between opioid use states. The definitions of "active" and "injection" opioid use can vary (but must be pre-specified) depending on the users' needs and available information. In the RESPOND Massachusetts base case, "active" opioid use is defined as any reported use in the previous seven days. "Injection" opioid use reflects any injection in the preceding seven days (a person who is both injecting and using oral opioids would be categorized as "injection" in RESPOND).
The care delivery module (Supplemental Figure 2) of RESPOND simulates OUD treatment and includes four treatment types: 1) outpatient buprenorphine (Bup), 2) outpatient injectable naltrexone (Ntx), 3) outpatient methadone (Mmt) maintenance, and 4) inpatient acute drug detoxification (detox). The model is adaptable to additional intervention types to better reflect local conditions and evolutions in the treatment field. For this study, a fifth block was added to capture the state of Corrections. In general, treatment episodes tend to decrease movement into active drug use, increase movement into non-active drug use, and have an independent effect on overdose rates conditional on active drug use. When population disengages from a treatment and is lost to follow-up, those people enter a corresponding "post-treatment state". The post-treatment state is a fixed interval during which relapse to active drug use is high, tolerance to opioids is lower than before treatment, and the risk of drug overdose among those actively using opioids is higher than it is in the no treatment state. The post-treatment state represents the period of vulnerability and excess overdose observed in real-world settings among patients who have recently relapsed to opioid use after a period of sustained abstinence.
The mortality module simulates both drug-related and competing risks deaths. RESPOND simulates overdose mortality by first simulating overdose incidence as a function of age and type of drug use (injection vs. non-injection use). Next, the model simulates a probability of death conditional on having had an opioid overdose. The model simulates competing causes of death through the use of standardized mortality ratios that are a function of age, sex, and type of opioid use (injection vs. non-injection).
The primary model outputs are: 1) All-cause mortality, 2) Overdose mortality, and 3) Number of people on treatment.
The simulation process is as follows: At simulation start, the model initiates a cohort of people currently living with OUD in the jurisdiction of interest. Based on data from that jurisdiction, the model assigns the current population to a drug use state, as well as a treatment block, such that the simulated population, including the prevalence of OUD treatment, reflects the status quo. Moving forward through simulated time, the sequence of simulation steps are: 1) aging of the population, 2) arrival of new population, 3) transition between OUD drug use states, 4) transitions into and out of treatment, 5) overdose, and 6) death.
At the end of this sequence of processes, the model advances simulated time by one cycle (week) and repeats the process. The simulation continues until a time horizon assigned by the user.
Full model schematics are available at https://www.syndemicslab.org/respond. Transitions between drug use compartments impact four important outcomes: 1) risk of overdose, 2) risk of death from competing causes, 3) health care utilization (cost), and 4) quality of life.
The primary sources of data for substance use transitions are studies from the medical literature.

Structural Assumptions:
• OUD is a remitting and relapsing process over a lifetime. There is no health state of OUD cure or permanent recovery.
• Transitions between OUD health states are not time updated.

Methodological Notes:
Supplemental Table 1 presents the key parameters related to OUD transitions for no treatment.
• CIs for proportions pE, pF, pB, pH, and pD are calculated using the normal approximation to binomial proportions.
• CIs for rates RA, RC, and RD are provided from the manuscript and calculated assuming Poisson distribution.
• Weekly rates and proportions, calculated from the respective overall estimates, are converted to weekly transition probabilities as indicated in the "Method" column in Supplemental Table 1

PG Probability of non-active injection to active injection
• p * C and p * G indicate the percentage of people relapsed within a month of discharge (after inpatient detoxification).
• The denominators for calculating weekly probabilities depend on whether the respective available proportion or rate estimates are yearly or monthly.

B.1.2 Overdose
Every person who is actively using opioids faces the risk of overdose. The probability of overdose depends on age, sex, and route of drug use (injection vs. non-injection). The simulation has no memory of past overdose events and does not include an elevated risk of repeat overdose after experiencing a first overdose event.
Structural Assumptions: • Experiencing overdose has no independent impact on current or future opioid use behaviors.
• Only the population that is in an active opioid use state faces the risk of overdose.
• The risk of overdose is different between no-treatment, treatment, and post-treatment episodes.
• The risk of overdose is lower while engaged in treatment compared to not-engaged, even among the population who are actively using drugs while engaged with treatment.
Methodological Notes: Counts of overdose are a target for model calibration. Supplemental Table 2  Note that RESPOND model simulations has weekly time cycles. Let (t) be the overdose rate in a specific block where = , or for blocks no-treatment, treatment, and post-treatment. Then, weekly overdose probabilities POD,B,t are calculated from the respective overdose rates as: where = , or .

B.1.2.2 No Treatment
Weekly overdose rates in no treatment (t) were calculated by applying a multiplier ∈ ℝ + on overall overdose rates ο(t) as (t) = ο(t) × . There were no data available to inform the rate multiplier . Therefore, we decided to calibrate overdose rate multiplier of no treatment block.

B.1.2.3 Overdose While on Treatment
The risk of overdose for people engaged in treatment, is derived by applying a multiplier parameter to the respective no-treatment ( ) estimates. i.e., the weekly overdose rate at time t for treatment is: where (t) is the no-treatment overdose rate at time t, and ∈ (0, 1).

B.1.2.4 Overdose During the Post-Treatment Period
During the post-treatment period, individuals face a risk of overdose higher than that of people who never initiated a treatment. Therefore, we model post-treatment overdose rates (t) with a multiplier (greater than 1) applied on no-treatment overdose rates such that (t) = (t) × . There were no data available to inform post-treatment overdose rate multiplier . Therefore, we decided to calibrate overdose rate multiplier of post-treatments.

B.1.3 Mortality
RESPOND simulates mortality through two independent mechanisms, fatal opioid overdose and nonoverdose death.

B.1.3.1 Fatal Overdose
The population that experiences overdose then faces a probability of death conditional on having had an overdose. This conditional probability of death, given an opioid overdose, is generalizable to all overdose cases and is therefore not stratified by age, sex, or OUD status. The population that survives an overdose does not change substance use as a result of the overdose. The probability of death conditional on having experienced an overdose is a time updated variable, reflecting changes to drug supply over time.
Adjusting the conditional probability of overdose death provides a mechanism to reflect the growing penetration of fentanyl in local drug supplies, which is a major dynamic underlying mounting overdose deaths in the U.S.
The probability ( ) of fatal overdose at year t is calculated as: where NFOD,t is the total number of fatal overdoses, and NOD,t is the total number of all-type overdoses.

B.1.3.2 Competing risks of death (non-overdose mortality)
Competing risks mortality includes deaths from conditions such as infectious endocarditis and sepsis, as well as medical comorbidities that accrue over a lifetime. The general approach to estimating competing risks of death is to apply standardized mortality ratios (SMRs) reflecting elevated mortality among drug users to age-sex stratified actuarial lifetables for the U.S. 19 .
where SMR is calculated as: where NNon-FOD is the number of observed deaths not due to opioid overdose, RD is the census death rate, and NOUD is the size of the OUD population based on chapter 55 estimation.
We construct CIs around the SMRs estimates assuming that the NNon-FOD follows a Poisson distribution and using the normal approximation (Supplemental Table 7 Finally, we convert weekly non-overdose death rates to probabilities as 1 − e − .

B.1.3.3 Summary of the combined impact of medications for opioid use disorder on all-cause mortality in RESPOND
Medications for opioid use disorder (MOUD) have two independent effects on mortality that combine to provide synergies in the simulation: 1. The population that is engaged with MOUD treatment experiences a net movement toward nonactive drug use. Because there is no risk of overdose while not using drugs, MOUD tend to decrease the rate of overdose in the population. In addition, movement out of active drug use states reduces exposure to the high standardized mortality ratios (SMRs) of active drug use and thereby reduce non-overdose mortality as well.
2. Among those who are actively using drugs when taking an MOUD, the MOUD has an independent effect on overdose risk, such that even those who are using have lower risk of death than those who are using drugs while not engaged with MOUD treatment.

B.1.4 Care Delivery
RESPOND models OUD while engaged with treatment using the same 4-state opioid use simulation that it uses to model OUD without treatment. The 4-state OUD simulation is embedded within all treatment episodes (blocks), such that individuals may both remain engaged with treatment, but also experience periods of drug use relapse. Each treatment type has its own bi-directional transition probabilities between active and non-active use. The net movement between active and non-active use while engaged with treatment favors movement to non-active use over time. The population that is lost to follow-up (disengages from care) must pass through a "post-treatment period" before rejoining the simulation of OUD. The post-treatment period is a four-week time, immediately following discontinuation of a treatment, during which the risk of relapse to drug use is high, as is the risk of overdose. Population that survives the post-treatment period transitions back to the simulation of OUD without treatment.

B.1.4.1 Movement From No-Treatment to Treatment Episodes
Structural Assumptions: • Only population in active opioid use states seeks OUD treatment. Population that is not currently using opioids does not seek treatment.
The main source of data to inform the probability of transition from no treatment to a treatment episode is the MA PHD.
Methodological Notes: Let denotes the weekly transition rates from no-treatment to treatment. Then, weekly transition probability from no-treatment to treatment is calculated from as: N Total,NoTrt → Trt .T : the total number of people with OUD "at risk" of transitioning from no-treatment to treatment T in January 2013 The weekly transition probability P NoTrt → Trt .T is estimated using data from the MA PHD repository 1 , and is stratified by age (16 groups: 5-year age-groups from 10-85, and >85years old), sex, and treatment (T= Detox, Mmt, Ntx, and Bup) (Supplemental Table 8

B.1.4.2 Treatment Initiation Effect
When population begins a treatment for opioid use disorder, for example out-patient buprenorphine, a portion of the population immediately transitions from active to non-active use. Following that initial "treatment initiation effect" there is bidirectional movement between active and nonactive use states, even while engaged with treatment. The main source of data for the treatment initiation effect and for substance use transitions while engaged with buprenorphine, naltrexone, or methadone is the NIDA CTN urine toxicology data. The CTN trials collected routine periodic urine toxicology from all participants.
While the published clinical trials results censored participants at the first relapse to drug use (the primary outcome of that trial), the trials continued to collect data from patients who experienced a relapse, such that the database includes longitudinal urine toxicology from patients who relapsed to active use, as well as some who remitted back to non-active use over the course of the trial. We analyzed those data in an "as treated" manner, such that RESPOND estimates realistic movements between active and non-active drug use states among people who are taking a medication. Note that relapsing to active drug use is not the same thing as loss to follow-up from treatment (see below).
Methodological Notes: Upon entering treatment, a proportion of the population immediately transitions from active to nonactive opioid use (Supplemental Table 9). This proportion pInit_Act→NonAct is stratified by treatment episode as follows: • methadone (Mmt): 0.57, based on the proportion of observed negative urine samples at week 5 We assume a binomial distribution and we use the Wald's method to calculate 95% CIs for the proportion pAct→NonAct representing the block initiation effect.
The model takes this input as the reciprocal, or pAct→Act. It is the direct reciprocal because the only transitions that can be made in this specific instance are Active to Non-active within the initial OUD type (injection vs. non-injection). There is no way to cross over to the other OUD type (ie there is no transition possible from Active non-injection to Non-active injection). These probabilities are listed below in Table   6. For example, take an individual who is in the 'Active non-injection' opioid use state and the 'No Treatment' care delivery state. When they move to the 'Buprenorphine' (Bup) care delivery state (based on their probability of movement identified in Table 5) they have a probability of 0.257 of immediately shifting to a 'Non-active non-injection' opioid use state. After that first week, their opioid use state transitions will instead be governed by the probabilities in Supplemental Table 10. *** There is no block initiation effect for population that is not currently using opioids, because only population that is currently using opioids seeks care in the model. **** Calibration parameter.

B.1.4.3 Transitions Between Active and Non-Active Opioid Use While Engaged with Treatment
We estimated Weekly OUD transition probabilities pTrt_Act→NonAct using Multi-State Models (MSMs) 15 . We fit separate models for each treatment: buprenorphine (Bup), naltrexone (Ntx), and methadone (Mmt), using data from the National Institute of Drug Abuse Clinical Trials Network (NIDA CTN) 4-6 .
Structural Assumptions: • Population engaged with treatment may move between active and non-active opioid use, but the population engaged with treatment does not change the route of administration of their opioid use. In other words, population that entered treatment using non-injection opioids will not escalate to injection drug use while still engaged with treatment (Core Simulation within OUD treatment episodes (blocks) -(eFigure ).
Transition probabilities between active and non-active states are the same for both injection and noninjection drug use. This structural assumption is confirmed by the MSM estimates for buprenorphine and methadone models, in which route was included as a model covariate, but was not a significant predictor of transition rates Methodological Notes: • Each MSM includes age and sex as covariates.
• Age is included as a continuous covariate in the MSM model, thus allowing estimation of the transition probabilities for age bins in which data are not available. We consider five 5 age groups: 10-19, 20-24, 25-34, 35-49, and 50-99 years old.
• OUD transition for Buprenorphine and Methadone: We keep all the weekly MSM estimates of OUD transition probabilities except week 1, which is considered as block initiation.
• OUD transition for Naltrexone: We delete the estimates for the first 4 weeks due to the inaccurate results from detoxification. Week 5 is also excluded from the analysis, as it is considered as block initiation.
• Transition probabilities from non-active to active use are defined as: * Probabilities listed here are transition probabilities to "Active" state. For ex:when initial OUD status is "Active", probabilities listed here are probabilities of staying in the "Active" state.

B.1.4.4 Probability of Loss to Follow-Up
In every time step, the population that is engaged with treatment faces a risk of disengaging from care and being lost to follow-up. Loss to follow-up differs from relapse to active drug use while remaining engaged with opioid treatment. The population that disengages with care and is lost to follow-up enters the "post-treatment state," during which time they have a high rate of relapse to active use and a high rate of overdose among active users. The post-treatment block represents the period of time immediately following discontinuation of a medication or release from an abstinence-based setting (acute drug detoxification center, residential drug treatment, or jail), when opioid tolerance is low and the risk of overdose is higher than that of a person who never initiated treatment.
The main source of data for estimating the probability of loss to follow-up is Market Scan, a large insurance claims database containing millions of individuals who have commercial insurance coverage.
As a randomized controlled trial, the CTN data cannot provide estimates of retention in care or loss to follow-up in the real world. We have previously published rates of loss to follow-up from buprenorphine and naltrexone treatment 16 . We therefore turn to Market Scan, which is nationally representative and reflects real-world practice in the U.S.
Structural Assumptions: • In RESPOND, the only way to transition into a post-treatment episode is from a corresponding treatment episode.
• The "No Treatment" block does not have a post-treatment episode.
• RESPOND also considers the probability of immediate relapse to active opioid use upon being lost to follow-up from treatment: Methodological Notes: The weekly transition probability from treatment to post-treatment is calculated as:  To simulate the demography and OUD epidemiology in the underlying population, RESPOND requires the initial cohort to be specified as follows: 1) age and sex distributions of people with opioid use disorder, 2) proportion of people beginning in each drug-use state, and 3) proportion of people within each treatment episode.
Structural Assumptions: • No population begins the simulation in a post-treatment block.
• RESPOND does not characterize the population by race or ethnicity.
Methodological Notes: Supplemental Table 12 presents the key parameters related to cohort initialization. The capture-recapture approach provides a method to estimate the total population with high-risk opioid use in a given calendar year, including those who have not been identified as being an opioid user and do not appear in medical claims or prevalence surveys. The previous work in provides estimates for year 2012, however, as it has been noted in the paper, the data sources used before 2013 were not complete 7 .
Hence, we used data in the years with more comprehensive data sources to predict the data in 2012 which has fewer data sources.
Second, we obtained the age and sex stratified counts for alive OUD population at the end of each year by subtracting death with opioid overdose involved from the estimated total size.  RESPOND simulates discrete time steps (rather than continuous time) and categorical age groups or "brackets" over the lifetime. The user can define the bounds of age groups so as to match the structure of the underlying population. Aging occurs as the population progresses to the next age group after a number of cycles that is determined by the size of the age brackets. The model employs a half-cycle correction and aging occurs in discrete steps, namely only at multiples of the age group size.
Structural Assumptions: • The entire population of the last age bracket (95 to 100-year-olds) is removed from the simulation at each aging cycle and replaced by the population from the previous age bracket. to the population occur at every time step. Arrival rates are stratified by age and sex and vary over time to realistically reflect the incidence of new opioid use disorder and/or movement into the jurisdiction.
Structural Assumptions: • All new populations enter the first block ("No Treatment" episode) and the first OUD state (currently active, non-injection).
• All new arrivals enter the simulation as active non-injectors under the "No Treatment" block but can transition to other feasible OUD states in the OUD transition module.
Methodological Notes: Supplemental Table 13 and Supplemental Table 14 present the key parameters related to entering cohort.
• The size of the entering cohort (Nenter,t) each year defined as: where NOUD,t is the estimated OUD population in year t, NFOD,t-1 is the estimated number of fatal overdoses, and NNon_FOD,t is the estimated total number of deaths from other causes, in the previous year (t-1).